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Abstract 


A  new  maximum  angle  method  has  been  developed  to  detennine  surface  mixed-layer  (a 
general  name  for  isothermal/constant-density  layer)  depth  from  profile  data.  It  has  three  steps: 
(1)  fitting  the  profile  data  with  a  first  vector  (pointing  downward)  from  a  depth  to  an  upper  level 
and  a  second  vector  (pointing  downward)  from  that  depth  to  a  deeper  level,  (2)  identifying  the 
angle  (varying  with  depth)  between  the  two  vectors,  (3)  finding  the  depth  (i.e.,  the  mixed  layer 
depth)  with  maximum  angle  between  the  two  vectors.  Temperature  and  potential  density  profiles 
collected  from  two  seagliders  in  the  Gulf  Stream  near  Florida  coast  during  14  November  -  5 
December  2007  were  used  to  demonstrate  its  capability.  The  quality  index  (1.0  for  perfect 
identification)  of  the  maximum  angle  method  is  about  0.96.  The  isothermal  layer  depth  is 
generally  larger  than  the  constant-density  layer  depth,  i.e.,  the  barrier  layer  occurs  during  the 
study  period.  Comparison  with  the  existing  difference,  gradient,  and  curvature  criteria  shows  the 
advantage  of  using  the  maximum  angle  method.  Uncertainty  due  to  varying  threshold  using  the 
difference  method  is  also  presented. 
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1.  Introduction 


Transfer  of  mass,  momentum,  and  energy  across  the  bases  of  surface  isothermal  and 
constant-density  layers  provides  the  source  for  almost  all  oceanic  motions.  Underneath  the 
surface  constant-density  and  isothermal  layers,  there  exist  layers  with  strong  vertical  gradients 
such  as  the  pycnocline  and  thermocline.  The  mixed  layer  depth  (MLD)  (a  general  name  for 
isothermal/constant-density  layer  depth)  is  an  important  parameter  which  influences  the 
evolution  of  the  sea  surface  temperature.  The  isothennal  layer  depth  (Hr)  is  not  necessarily 
identical  to  the  constant-density  layer  depth  (Hd)  due  to  salinity  stratification.  There  are  areas  of 
the  World  Ocean  where  HT  is  deeper  than  HD  (Lindstrom  et  ah,  1987;  Chu  et  ah,  2002;  de  Boyer 
Montegut  et  ah,  2007).  The  layer  difference  between  Hd  and  Hr  is  defined  as  the  barrier  layer. 
For  example,  barrier  layer  was  observed  from  a  seaglider  in  the  western  North  Atlantic  Ocean 
near  the  Florida  coast  (30.236°N,  78.092°W)  at  GMT  23:20  on  November  19,  2007  (Fig.  1).  The 
barrier  layer  thickness  (BLT)  is  often  referred  to  as  the  difference,  BLT  =  Hr  -  Hd .  Less 
turbulence  in  the  barrier  layer  than  in  the  constant-density  layer  due  to  strong  salinity 
stratification  isolates  the  constant-density  layer  water  from  cool  thennocline  water.  This  process 
regulates  the  ocean  heat  budget  and  the  heat  exchange  with  the  atmosphere,  and  in  turn  affects 
the  climate. 

Objective  and  accurate  identification  of  Hr  and  Hd  is  important  for  the  determination  of 
barrier  layer  occurrence  and  its  climate  impact.  Three  types  of  criteria  (difference,  gradient,  and 
curvature)  are  available  for  identifying  I  It  from  profiling  data.  The  difference  criterion  requires 
the  deviation  of  T  (or  p)  from  its  near  surface  (i.e.,  reference  level)  value  to  be  smaller  than  a 
certain  fixed  value.  The  gradient  criterion  requires  dT/  dz  (or  dp/  8 z)  to  be  smaller  than  a  certain 
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fixed  value.  The  curvature  criterion  requires  82T/  dz2  (or  d2p/dz2)  to  be  maximum  at  the  base  of 
mixed  layer  (z  =  -Ho). 

The  difference  and  gradient  criteria  are  subjective.  For  example,  the  criterion  for 
detennining  HT  for  temperature  varies  from  0.8°C  (Kara  et  ah,  2000),  0.5°C  (Wyrtki,  1964)  to 
0.2°C  (de  Boyer  Montegut  et  ah,  2004).  The  reference  level  changes  from  near  surface  (Wyrtki, 
1964)  to  10  m  depth  (de  Boyer  Montegut  et  ah,  2004).  The  criterion  for  determining  Ho  for 
potential  density  varies  from  0.03  kg/m  (de  Boyer  Montegut  et  ah,  2004),  0.05  kg/m  (Brainerd 
and  Gregg,  1995),  to  0.125  kg/m3  (Suga  et  al.,  2004).  Defant  (1961)  was  among  the  first  to  use 
the  gradient  method.  He  used  a  gradient  of  0.015°C/m  to  determine  / //  for  temperature  of  the 
Atlantic  Ocean;  while  Lukas  and  Lindstrom  (1991)  used  0.025°C/m.  The  curvature  criterion  is 
an  objective  method  (Chu  et  al,  1997,  1999,  2000;  Lorbacher  et  al.,  2006);  but  is  hard  to  use  for 
profile  data  with  noise  (even  small),  which  will  be  explained  in  Section  5.  Thus,  it  is  urgent  to 
develop  a  simple  objective  method  for  determining  MLD  with  the  capability  of  handling  noisy 
data. 

In  this  study,  a  new  maximum  angle  method  has  been  developed  to  determine  ///  and  Ho 
and  the  gradients  of  the  thennocline  and  pycnocline  from  profiles  and  tested  using  data  collected 
by  two  seagliders  of  the  Naval  Oceanographic  Office  in  the  Gulf  Stream  region  near  the  Florida 
coast  during  14  November  -  5  December  2007,  with  comparison  to  the  existing  methods.  The 
results  demonstrate  its  capability.  The  outline  of  this  paper  is  listed  as  follows.  Section  2 
describes  hydrographic  data  from  the  two  seagliders.  Section  3  presents  the  methodology. 
Sections  4  and  5  compare  the  maximum  angle  method  with  the  existing  methods.  Section  6 
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shows  the  occurrence  of  barrier  layer  in  the  western  North  Atlantic  Ocean.  Section  7  presents 
the  conclusions. 

2.  Seaglider  Data 

Two  seagliders  were  deployed  in  the  Gulf  Stream  region  near  the  Florida  coast  by  the 
Naval  Oceanographic  Office  (Mahoney  et  ah,  2009)  from  two  nearby  locations  on  14  November 
2007  with  one  at  (79.0°  W,  29.5°N)  (Seaglider-A)  and  the  other  at  (79.0°  W,  29.6°N)  (Seaglider- 
B)  (Fig.  2).  Seaglider-A  (solid  curve)  moved  toward  northeast  to  (78.1°W,  30.25°N),  turned 
anticyclonically  toward  south  and  finally  turned  cyclonically  at  (78.4°W,  29.6°N).  Seaglider-B 
(dashed  curve)  moved  toward  north  to  (79.0°W,  30.0°N),  turned  toward  northeast  and  then 
anticyclonically,  and  finally  turned  cyclonically. 

The  seaglider  goes  up  and  down  in  oblique  direction,  not  vertical.  Data  collected  during  a 
downward-upward  cycle  are  divided  into  two  parts  with  the  first  one  from  the  surface  to  the 
deepest  level  and  the  second  one  from  the  deepest  level  to  the  surface.  Each  part  represents  an 
individual  profile  with  the  averaged  longitude  and  latitude  of  the  data  points  as  the  horizontal 
location.  Such  created  temperature  and  potential  density  profile  data  went  through  quality 
control  procedures;  such  as,  min-max  check  (e.g.,  disregarding  any  temperature  data  less  than  - 
2°C  and  greater  than  35°C),  error  anomaly  check  (e.g.,  rejecting  temperature  and  salinity  data 
deviating  more  than  3  standard  deviations  from  climatology),  seaglider-tracking  algorithm 
(screening  out  data  with  obvious  seaglider  position  errors),  max-number  limit  (limiting  a 
maximum  number  of  observations  within  a  specified  and  rarely  exceeded  space-time  window), 
and  buddy  check  (tossing  out  contradicting  data).  The  climatological  data  set  used  for  the  quality 
control  is  the  Navy’s  Generalized  Digital  Environmental  Model  (GDEM)  climatological 
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temperature  and  salinity  data  set.  After  the  quality  control,  514  profiles  of  (T,  p )  are  available 
with  265  profiles  from  Seaglider-A  and  249  profiles  from  Seaglider-B.  The  vertical  resolution  of 
the  profile  varies  from  less  than  1  m  for  upper  10  m,  to  approximately  1-3  m  below  10  m  depth. 
All  the  profiles  are  deeper  than  700  m  and  clearly  show  the  existence  of  layered  structure: 
isothennal  (constant-density)  layer,  thermocline  (pycnocline),  and  deep  layer. 

3.  Determination  of  ( HT ,  HD) 

Let  potential  density  and  temperature  profiles  be  represented  by  [p(zk),  T(zk)].  Here,  k 
increases  downward  with  k  =  1  at  the  surface  or  nearest  to  the  surface  and  K  the  total  number  of 
the  data  points  for  the  profile.  The  potential  density  profile  is  taken  for  illustration  of  the  new 
methodology.  Let  (pmax  ,  Pmm)  be  the  maximum  and  minimum  of  the  profile  p(zk).  Starting  from 
zi  downward,  depth  with  pmm  (zmm)  and  depth  with  pmax  ( zmax )  are  found.  Without  noise,  zmm 
should  be  z\  and  zmax  should  be  zK.  The  vertical  density  difference,  A p  =  pmax  -  pmin,  represents 
the  total  variability  of  potential  density.  Theoretically,  the  variability  is  0  in  the  constant-density 
layer  and  contains  large  portion  in  the  pycnocline.  It  is  reasonable  to  identify  the  main  part  of  the 
pycnocline  between  the  two  depths:  Z(0.p  with  0.1  Ap  and  z((k7)  with  of  0.7A p  relative  to  pmm  (Fig. 
3).  Let  n  be  the  number  of  the  data  points  between  z(o.i)  and  z(  0.7)  and  m  =  min(n,20)  . 

At  depth  Zk  (marked  by  a  circle  in  Fig.  4),  a  first  vector  (Ai,  downward  positive)  is 
constructed  with  linear  polynomial  fitting  of  the  profile  data  from  z/<:/  to  zk  with 

k- 1,  for  k  <  m 

j  =  {  ’  7  (1) 

m,  tor  k  >  m 

A  second  vector  (A2,  pointing  downward)  from  one  point  below  that  depth  (i.e.,  zk+\)  is 

constructed  to  a  deeper  zk+m.  The  dual-  linear  fitting  can  be  represented  by 
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(2) 


P(z)  =  { 


Zk-j’Zk-j+l’—Zk 

Zk+l’—Zk+m 


where  ck\  c(k2>,  Gk\  G[2)  are  the  fitting  coefficients. 


At  the  constant-density  (isothermal)  layer  depth,  the  angle  Ok  reaches  its  maximum  value 
if  Zk  is  the  MLD  (see  Fig.  4a),  and  smaller  if  z/c  is  inside  (Fig.  4b)  or  outside  (Fig.  4c)  of  the 
mixed  layer.  Thus,  the  maximum  angle  principle  can  be  used  to  detennine  the  mixed  (or 
isothermal)  layer  depth, 

0k  — >  max,  HD=-zk. 


In  practical,  the  angle  0 k  is  hard  to  calculate,  so  tan  Ok  is  used  instead,  i.e., 


tan  0k  — »  max,  II n 


—z 


G0>  =  G)'\  GK2)=G\ 


;(l) 


k ’ 


(2) 


;(2) 


(3) 


where  G(1)  «  0  ,  is  the  vertical  gradient  in  the  mixed  layer;  and  G(2)  is  the  vertical  gradient  in  the 
thermocline  (pycnocline).  With  the  given  fitting  coefficients Gk\  G{k\  tan  Ok  can  be  easily 
calculated  by 


_  G^-G? 
l  +  G^G™  ‘ 


(4) 


The  maximum  angle  method  [i.e.,  (1)-  (4)]  was  used  to  calculate  HT  and  HD  from  514  pairs  of 
temperature  and  potential  density  profiles  from  the  two  seagliders.  With  high  vertical  resolution 
of  the  data,  Hr  and  Hd  can  be  detennined  for  all  profiles.  The  potential  density  profile  at  the 
station  (shown  in  Fig.  1)  located  at  (30.236°N,  78.092°W)  is  taken  as  an  example  for  illustration 
(Fig.  5a).  At  z  =  -HD,  tan  Ok  has  maximum  value  (Fig.  5b). 
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Advantage  of  the  maximum  angle  method  is  described  as  follows.  Different  from  the 
existing  methods,  the  maximum  angle  method  not  only  uses  the  main  feature  (vertically  uniform) 
of  mixed  layer  such  as  in  the  difference  and  gradient  criteria,  but  also  uses  the  main 
characteristics  (sharp  gradient)  below  the  mixed  layer  (see  Fig.  3).  After  MLD  is  determined,  the 
vertical  gradient  of  the  thermocline  (pycnocline),  G(2),  is  also  calculated.  The  dataset  of  G(2)  is 
useful  for  studying  the  heat  and  salt  exchange  between  the  ocean  upper  and  lower  layers. 
Besides,  the  maximum  angle  method  is  less  subjective  comparing  to  the  existing  methods.  The 
only  external  parameters  (10%  and  70%)  are  used  to  detennine  zq.\  and  zoj,  and  in  turn  to 
determine  the  length  of  the  vectors  Ai  and  A2. 

Disadvantage  of  the  maximum  angle  method  is  due  to  the  use  of  two  linear  regressions 
[see  Eq(2)].  Reliable  regression  needs  sufficient  sample  size.  For  profiles  with  very  few  data 
points  (low  resolution),  the  maximum  angle  method  might  not  work.  The  seaglider  data 
described  in  Section  2  are  high-resolution  profiles,  and  therefore  are  perfect  for  the  test  of  the 
maximum  angle  method. 

4.  Comparison  between  Maximum  Angle  and  Threshold  Methods 

Lorbacher  et  al.  (2006)  proposed  a  quality  index  (Ql) 

rmsd(A-A)U^)  m 

rmsd(/9Jt  ~ Pk)  \(h,,\.5xhd)  ’ 

to  evaluate  various  schemes  for  HD  (or  ///)  determination.  Here,  pk  =  p(zk),  is  the  observed 
profile,  pk  =mean  potential  density  between  z\  and  z  =  -HD.  For  a  perfect  identification, 

rmsd ( pk  ~ pk)  \(Hi  hd)=  0,  QI  =  1.  The  higher  the  Ql,  the  more  reliable  identification  of  MLD 


8 


would  be.  Usually,  HD  is  defined  with  certainty  if  QI  >  0.9;  can  be  determined  with  uncertainty  if 
0.9  >  QI  >  0.5;  and  can’t  be  identified  if  QI  <  0.5. 

Capability  of  the  maximum  angle  method  is  demonstrated  through  comparison  with  the 
existing  threshold  method.  Since  the  MLD  based  on  a  difference  criterion  is  more  stable  than  the 
MLD  based  on  a  gradient  criterion,  which  requires  sharp  gradient-resolved  profiles  (Brainerd 
and  Gregg,  1995),  the  difference  threshold  method  is  used  for  the  comparison.  Four  sets  of 
isothennal  depth  were  obtained  from  514  temperature  profiles  of  the  two  seagilders  using  the 
maximum  angle  method,  0.2°C  (de  Boyer  Montegut  et  al.,  2004),  0.5°C  (Monterey  and  Levitus, 
1997),  and  0.8°C  (Kara  et  ah,  2000)  thresholds.  Fig.  6  shows  the  histograms  of  514  Hr  -values 
for  the  four  methods.  Difference  of  the  histograms  among  0.2°C  (Fig.  6b),  0.5°C  (Fig.  6c),  and 
0.8°C  (Fig.  6d)  thresholds  implies  uncertainty  using  the  difference  method.  Table- 1  shows  the 
statistical  characteristics  of  Hr  determined  by  the  four  methods.  The  mean  (median)  Hr  value  is 
77.2  m  (73.2  m)  using  the  maximum  angle  method.  For  the  difference  method,  it  increases  with 
the  value  of  the  threshold  from  71.9  m  (71.8  m)  using  0.2°C,  82.0  m  (77.9  m)  using  0.5°C,  to 
87.6  m  (82.9  m)  using  0.8°C. 

The  Gaussian  distribution  has  skewness  of  0  and  kurtosis  of  3.  Obviously,  the  four 
histograms  show  non-Gaussian  features.  ///  is  positively  skewed  when  it  is  identified  using  all 
the  four  methods.  The  skewness  of  Hj  is  sensitive  to  the  thresholds:  0.21  using  0.2°C,  1.13 

using  0.5°C,  and  1.25  using  0.8°C.  It  is  0.69  using  the  maximum  angle  method.  The  kurtosis  of 
Ht  is  larger  than  3  for  all  the  four  methods  and  sensitive  to  the  thresholds:  3.48  using  0.2°C,  4.46 
using  0.5°C,  and  4.35  using  0.8°C.  It  is  3.80  using  the  maximum  angle  method. 
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The  histograms  of  514  QI -values  are  negatively  skewed  for  the  four  methods  (Fig.  7). 
Most  Ql-values  are  larger  than  0.980  with  a  mean  value  of  0.965  using  the  maximum  angle 
method  (Fig.  7a)  and  are  lower  using  the  threshold  method  (Figs.  6b-d).  The  mean  Ql-value 
reduces  from  0.881  with  0.2°C  threshold  (Fig.  7b),  0.858  with  0.5°C  threshold  (Fig.  7c),  to  0.833 
with  0.8°C  threshold  (Fig.  7d). 

Uncertainty  of  the  difference  method  from  one  to  another  threshold  can  be  identified  by 
computing  the  relative  root-mean  square  difference  (RRMSD), 

i  n  * 

RRMSD  =  ,  (6) 

where  N  =  5 14,  is  the  number  of  the  seaglider  profiles;  ( H\X) ,  H)2) ,  i  =  1,  2,  . . .,  N)  are  two  sets 

of  MLD  identified  by  the  difference  method  using  two  different  criteria.  The  RRMSD  of  Hr  is 
20.1%  between  0.2°C  and  0.5°C  thresholds,  28.5%  between  0.2°C  and  0.8°C  thresholds,  and 
10.0%  between  0.5°C  and  0.8°C  thresholds. 

Similarly,  four  sets  of  constant-density  depth  were  obtained  from  514  potential  density 
profiles  of  the  two  seagilders  using  the  maximum  angle  method,  0.03  kg/m  (de  Boyer  Montegut 
et  ah,  2004),  0.05  kg/m3  (Brainerd  and  Gregg,  1995),  and  0.125  kg/m3  (Monterey  and  Levitus, 
1997)  thresholds.  Fig.  8  shows  the  histograms  of  514  Hd  -values  for  the  four  methods.  I  Id  is 
positively  skewed  when  it  is  identified  using  the  maximum  angle  method.  Difference  of  the 
histograms  among  0.03  kg/m3  (Fig.  8b),  0.05  kg/m3  (Fig.  8c),  to  0.125  kg/nr  (Fig.  8d)  threshold 
implies  uncertainty  in  the  difference  method.  Table-2  shows  the  statistical  characteristics  of  Hd 
determined  by  the  four  methods.  The  mean  (median)  Hd  value  is  73.2  m  (70.4  m)  using  the 
maximum  angle  method.  It  increases  with  the  value  of  the  threshold  from  53.3  m  (60.9  m)  using 
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0.03  kg/m3,  59.3  m  (66.2  m)  using  0.05  kg/m3,  to  68.0  m  (71.6  m)  using  0.125  kg/m3.  The 
skewness  of  Hd  is  slightly  positive  (0.28)  when  it  is  identified  using  the  maximum  angle  method 
and  slightly  negative  when  it  is  identified  using  the  threshold  method.  The  negative  skewness 
enhances  with  the  threshold  from  -0.06  using  0.03  kg/m  ,  -0.24  using  0.05  kg/m ,  to  -0.59 
using  0.125  kg/m3.  The  kurtosis  of  HD  is  4.32  using  the  maximum  angle  method,  and  varies  with 
the  threshold  when  the  difference  method  is  used.  It  is  2.37  using  0.03  kg/m  ,  2.74  using  0.05 
kg/m3,  and  3.67  using  0.125  kg/m3. 

The  histograms  of  514  QI -values  are  negatively  skewed  for  the  four  methods  (Fig.  9). 
Most  Ql-values  are  larger  than  0.980  with  a  mean  value  of  0.966  for  the  maximum  angle  method 
(Fig.  9a)  and  are  lower  using  the  threshold  method  (Figs.  8b-d)  comparing  to  the  maximum  angle 
method.  The  mean  Ql-value  reduces  from  0.837  with  0.03  kg/m3  threshold  (Fig.  9b),  0.859  with 
0.05  kg/m3  threshold  (Fig.  9c),  and  0.872  with  0.125  kg/m3  threshold  (Fig.  9d).  The  RRMSD 
of  Hd  is  29.3%  between  0.03  kg/m  '  and  0.05  kg/m3  thresholds,  44.7%  m  between  0.03  kg/m3  and 
0. 125  kg/m3  thresholds,  and  27.8%  between  0.05  kg/m3  and  0. 125  kg/m3  thresholds. 

5.  Comparison  between  Maximum  Angle  and  Curvature  Methods 

Both  curvature  and  maximum  angle  methods  are  objective.  To  illustrate  the  superiority  of 
the  maximum  angle  method,  an  analytical  temperature  profile  with  / //  of  20  m  is  constructed 


T(z)  = 


21°C, 

2 1°  C  +  0.25°  Cx(z  + 20  m) 
6  z  +  40  nC 


7°C  +  9°C  x  exp 


50  m 


) 


-20  m  <  z  <  0  m 
-40  m  <  z  <  -20  m 


-100  m  <  z  <  -40  m. 


(V) 
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This  profile  was  discretized  with  vertical  resolution  of  1  m  from  the  surface  to  10  m  depth  and  of 
5  m  below  10  m  depth.  The  discrete  profile  was  smoothed  by  5-point  moving  average  to  remove 
the  sharp  change  of  the  gradient  at  20  m  and  40  m  depths.  The  smoothed  profile  data  \T(zk)\  is 
shown  in  Fig.  10a. 

The  second-order  derivatives  of  T(z/ :)  versus  depth  is  computed  by  nonhomogeneous 
mesh  difference  scheme, 


82T 

i 

(  T  —  T 

1k+ 1  1k 

Tk 

7 

i 

dz2 

z‘  ArVi 

V  Zk+ 1  _  Zk 

N 

~Zk-x) 

Here,  k  =  1  refers  to  the  surface,  with  increasing  values  indicating  downward  extension  of  the 
measurement.  Eq.(8)  shows  that  we  need  two  neighboring  values,  Tk.\  and  Tk+l,  to  compute  the 
second-order  derivative  at  zk .  For  the  surface  and  100  m  depth,  we  use  the  next  point  value,  that 
is, 


82T 

\  5271 1 

a2r. 

82T , 

dz2 

|  z=0  ^^2  z=~^  m , 

dz2  |z=- 

-mm-  dz2  \ 

Fig.  10b  shows  the  calculated  second-order  derivatives  from  the  profile  data  listed  in  Table  1. 

Similarly,  tan  9k  is  calculated  using  Eq.(4)  for  the  same  data  profile  (Fig.  10c).  For  the  profile 

2  2 

data  without  noise,  both  the  curvature  method  (i.e.,  depth  with  minimum  8  T/dz  ,  see  Fig.  10b) 
and  the  maximum  angle  method  [i.e.,  depth  with  max  (tan  9),  see  Fig.  10c)]  identified  the 
isothennal  depth,  i.e.,  Ht=  20  m. 

One  thousand  'contaminated'  temperature  profiles  are  generated  by  adding  random  noise 
with  mean  of  zero  and  standard  deviation  of  0.02°C  (generated  by  MATLAB)  to  the  original 
profile  data  at  each  depth.  An  example  profile  is  shown  in  Fig.  11a,  as  well  as  the  second  order 
derivative  (82T/8z2)  and  tan  9.  Since  the  random  error  is  so  small  (zero  mean,  0.02°C  standard 
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deviation,  within  the  instrument’s  accuracy),  we  may  not  detect  the  difference  between  Fig.  10a 
and  Fig.  1  la  by  eyes.  However,  the  isothermal  depth  is  9  m  (error  of  1 1  m)  using  the  curvature 
method  (Fig.  lib)  and  20  m  (Fig.  lid)  using  the  maximum  angle  method. 

Usually,  the  curvature  method  requires  smoothing  for  noisy  data  (Chu,  1999;  Lorbacher 
et  al.,  2006).  To  evaluate  the  usefulness  of  smoothing,  a  5-point  moving  average  was  applied  to 
the  1000  “contaminated”  profile  data.  For  the  profiles  (Fig.  11a)  after  smoothing,  the  second 
derivatives  were  again  calculated  for  the  smoothed  profiles.  The  isothermal  depth  identified  for 
the  smoothed  example  profile  is  8  m  (Fig.  1  lc).  Performance  for  the  curvature  method  (with 
and  without  smoothing)  and  the  maximum  angle  method  is  evaluated  with  the  relative  root-mean 
square  error  (RRMSE), 

r®mse  =  -T7JT|;(//?)-//;')!  ,  (io) 

where  //"'  (=  20  m)  is  for  the  original  temperature  profile  (Fig.  10a);  N  (=  1000)  is  the  number 

of  “contaminated”  profiles;  and  Hj]  is  the  calculated  for  the  /-th  profile.  Table  3  shows  the 

frequency  distributions  and  RRMSEs  of  the  calculated  isothennal  depths  from  the 
“contaminated”  profile  data  using  the  curvature  method  (without  and  with  5  point-moving 
average)  and  the  maximum  angle  method  without  smoothing.  Without  5-point  moving  average, 
the  curvature  method  identified  only  6  profiles  (out  of  1000  profiles)  with  Hr  of  20  m,  and  the 
rest  profiles  with  HT  ranging  relatively  evenly  from  1  m  to  10  m.  The  RRMSE  is  76%.  With  5- 
point  moving  average,  the  curvature  method  identified  413  profiles  with  HT  of  20  m,  164 
profiles  with  Hr  of  15  m,  3  profiles  with  Hr  of  10  m,  and  the  rest  profiles  with  ///  ranging 
relatively  evenly  from  2  m  to  8  m.  The  RRMSE  is  50%.  However,  without  5-point  moving 
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average,  the  maximum  angle  method  identified  987  profiles  with  Ht  of  20  m,  and  13  profiles 
with  Ht  of  15  m.  The  RRMSE  is  less  than  3%. 

6.  Existence  of  Barrier  Layer 

With  Hd  and  ///,  the  BLT  is  easily  calculated  from  all  5 14  profiles.  BLT  is  plotted  versus 
time  in  Fig.  12a  (Seaglider-A)  and  Fig.  12b  (Seaglider-B)  using  the  maximum  angle  method. 
These  two  figures  show  a  rather  frequent  occurrence  of  barrier  layer  in  the  western  North 
Atlantic.  For  example,  among  265  stations  from  Seaglider-A,  there  are  176  stations  where 
barrier  layer  occurs.  The  barrier  layer  occur  in  66.4%.  The  BLT  has  a  maximum  value  of  30.0  m 
on  30  November  2007.  Among  249  stations  from  Seaglider-B,  there  are  131  stations  where 
barrier  layer  occurs.  The  barrier  layer  occurs  in  52.6%.  In  this  l°x  2°  area,  variation  of  the  BLT 
has  complex  pattern  and  intennittent  characteristics. 

From  Tables  1  and  2,  the  mean  values  of  Ht  and  Hd  are  77.2  m  and  73.2  m,  which  lead  to 
the  mean  BLT  of  4.0  m.  When  the  difference  method  is  used,  identification  of  BLT  depends  on 
the  threshold.  For  example,  de  Boyer  Montegut  et  al.  (2004)  used  0.2°C  and  0.03  kg/m3.  From 
Tables  1  and  2,  the  mean  values  of  Ht  and  Hd  are  71.9  m  and  53.3  m,  which  lead  to  the  mean 
BLT  of  18.6  m.  Monterey  and  Levitus  (1997)  used  0.5°C  and  0.125  kg/m3.  From  Tables  1  and 
2,  the  mean  values  of  Ht  and  Hd  are  82.0  m  and  68.0  m,  which  lead  to  the  mean  BLT  of  14.0  m. 
Comparing  the  existing  difference  methods,  the  barrier  layer  has  less  chance  to  occur  using  the 
maximum  angle  method. 

7.  Conclusions 

A  new  maximum  angle  method  is  proposed  in  this  study  to  identify  isothermal  and 

constant-density  layer  depths.  First,  two  vectors  (both  pointing  downward)  are  obtained  using 
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linear  fitting.  Then,  the  tangent  of  the  angle  (tan  0)  between  the  two  vectors  is  calculated  for  all 
depth  levels.  Next,  the  isothermal  (or  constant-density)  depth  which  corresponds  to  the 
maximum  value  of  (tan  9)  is  found.  Two  features  make  this  method  attractive:  (a)  less  subjective 
and  (b)  capability  to  process  noisy  data.  The  temperature  and  potential  density  profiles  collected 
from  two  seagliders  in  the  Gulf  Stream  near  Florida  coast  during  14  November  -  5  December 
2007  were  used  for  evaluating  the  new  algorithm.  With  high  quality  indices  (QI  ~  96%),  the 
maximum  angle  method  not  only  identify  Hd  and  ///,  but  also  the  potential  density 
(temperature)  gradient  [G(2)]  below  z  =  -Hd  (z  =  -  ///).  Weakness  of  the  maximum  angle  method 
is  due  to  the  sample  size  requirement  of  the  regression.  For  low  resolution  profiles,  the  maximum 
angle  method  might  not  be  suitable. 

Uncertainty  in  detennination  of  (//?,  I  In)  due  to  different  thresholds  is  demonstrated 
using  the  same  seaglider  data.  Histogram  of  Hr  (Hd)  changes  evidently  when  different 
thresholds  are  used:  0.2°C  (0.03  kg/m3),  0.5°C  (0.05  kg/m3),  and  0.8°C  (0.125  kg/m3).  The 
RRMSD  of  Ht  is  20.1%  between  0.2°C  and  0.5°C  thresholds,  28.5%  between  0.2°C  and  0.8°C 
thresholds,  and  10.0%  between  0.5°C  and  0.8°C  thresholds.  The  RRMSD  of  HD  is  29.3% 
between  0.03  kg/m3  and  0.05  kg/m3  thresholds,  44.7%  m  between  0.03  kg/nr3  and  0.125  kg/m3 
thresholds,  and  27.8%  between  0.05  kg/nr  and  0.125  kg/m  thresholds.  Such  large  values  of 
RRMSD  make  the  difference  method  unreliable. 
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Table  1.  Statistical  characteristics  of  Ht  identified  from  the  two  seagliders  using  the  maximum 
angle,  0.2°C,  0.5°C,  and  0.8°C  thresholds. 


Maximum 

Angle 

0.2°C 

threshold 

0.5°C 

threshold 

0.8°C 

threshold 

Mean  (m) 

77.2 

71.9 

82.0 

87.6 

Median  (m) 

73.2 

71.8 

77.9 

82.9 

Standard 
Deviation  (m) 

18.3 

23.4 

18.4 

18.0 

Skewness 

0.69 

0.21 

1.13 

1.25 

Kurtosis 

3.80 

3.48 

4.46 

4.35 

Table  2.  Statistical  characteristics  of  Hd  identified  from  the  two  seagliders  using  the  maximum 
angle,  0.03  kg/nr,  0.05  kg/m3,  and  0.125  kg/m3  thresholds. 


Maximum 

Angle 

0.03  kg/m3 
threshold 

0.05  kg/m3 
threshold 

0.125  kg/m3 
threshold 

Mean  (m) 

73.2 

53.3 

59.3 

68.0 

Median  (m) 

70.4 

60.9 

66.2 

71.6 

Standard 
Deviation  (m) 

19.1 

32.9 

31.0 

28.4 

Skewness 

0.28 

-0.06 

-0.24 

-0.59 

Kurtosis 

4.32 

2.37 

2.74 

3.67 
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Table  3.  Frequency  distributions  and  RRMSE  of  calculated  isothennal  depths  from  the  data 
consisting  of  the  profile  (indicated  in  Table  1)  and  random  noise  with  mean  of  0  and  standard 
deviation  of  0.02°  C  using  the  curvature  method  (without  and  with  5  point-moving  average) 
and  the  maximum  angle  method  without  smoothing.  The  total  contaminated  data  profiles  are 


1000. 


Isothermal 

Layer 

Depth 

(m) 

Frequency  : 
Curvature 
(without 
smoothing) 

Frequency: 
Curvature 
(with  smoothing) 

Frequency 

Maximum  Angle 
(without  smoothing) 

1 

103 

0 

0 

2 

125 

83 

0 

3 

103 

55 

0 

4 

126 

44 

0 

5 

98 

52 

0 

6 

95 

47 

0 

7 

121 

43 

0 

8 

105 

96 

0 

9 

118 

0 

0 

10 

0 

3 

0 

15 

0 

164 

13 

20 

6 

413 

987 

Total 

1000 

1000 

1000 

RRMSE 

76% 

50% 

<3% 
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2007-11-19  23:19:45  lat:  30.236  Ion: -78.092 


Fig.  1.  Isothermal,  constant-density,  and  barrier  layers  were  observed  by  a  seaglider  in  the 
western  North  Atlantic  Ocean  near  the  Florida  coast  (30.236°N,  78.092°W)  at  GMT  23:20  on 
November  19,  2007. 
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Latitude 


Fig.  2.  (a)  Location  of  the  glider  data,  and  (b)  drifting  paths  of  two  gliders  with  the  marked 
station  by  circle  for  demonstration. 


(a)  (b) 


2007-11-19  23:19:45  lat:  30.236  Ion: -78.092 


density  (kg/m3)  tan(a) 


Fig.  5.  Determination  of  Ho  using  the  maximum  angle  method:  (a)  density  profile  at  the  seaglider 
station  (shown  in  Fig.  1)  located  at  (30.236°N,  78.092°W)  at  GMT  23:20  on  November  19,  2007, 
and  (b)  calculated  tan  oca.  It  is  noted  that  the  depth  of  the  maximum  tan  oca  corresponds  to  the 
mixed  layer  depth  and  only  the  upper  part  of  the  potential  density  profile  is  shown  here. 
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Frequency  Frequency 


0 


(a)  Analytic  Temperature  Profile 


(b)  Curvature  Method 


(c)  M-A  Method 


Fig.  10.  (a)  Smoothed  analytic  temperature  profile  (6)  by  5-point  moving  average,  calculated  (b) 
(< d2T/dz2)k ,  and  (c)  (tan  Of  from  the  profile  data  (Fig.  7a).  At  20  m  depth,  ( d2T/dz2)k  has  a 
minimum  value,  and  (tan  0)/(  has  a  maximum  value. 
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(a)  With  Noise  (a=0.01  C)  (b)  Curvature  Method  (c)  Smoothed  Curvature  Method  (d)  Maximum  Angle  Method 


Fig.  11.  One  out  of  1000  realizations:  (a)  temperature  profile  shown  in  Fig.  7a  contaminated  by 
random  noise  with  mean  of  zero  and  standard  deviation  of  0.02°C,  (b)  calculated  (d2T/dz2)ii  from 
the  profile  data  (Fig.  8a)  without  smoothing,  (c)  calculated  (3  773z  )*  from  the  smoothed  profile 
data  (Fig.  8a)  with  5-point  moving  average,  and  (d)  calculated  (tan  0)*  from  the  profile  data 
(Fig.  8a)  without  smoothing. 
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